Bending Wavelet for Flexural Impulse Response 



Richard BiissowO 
Einsteinufer 25, 10587 Berlin 

Richard BussovvQ 



o 
o 

(N 



(Dated: February 2, 2008) 

The work addresses the definition of a wavelet that is adapted to analyse a flexural impulse 
response. The wavelet gives the opportunity to directly analyse the dispersion characteristics 
of a pulse. The aim is to localize a source or to measure material parameters. An overview of 
the mathematical properties of the wavelet is presented. An algorithm to extract the dispersion 
characteristics with the use of genetic algorithms is outlined. The application of the wavelet is 
shown in an example and experiment. 
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I. INTRODUCTION 

The Morlet wavelet transform is a popular method for 
time-frequency analysis. Its application for acoustic sig- 
nals can be found in several publications. These publi- 
cations deal for example with the analysis of dispersive 
wavesi^, source or damage localizationMiiL&L ^ investi- 
gation of system parameter s] 9 ' 10 or active control^. A 
comparison of the short time Fourier transform and the 
Morlet wavelet transform is done by Kim et.al^. It is 
found that the continuous wavelet transform CWT of 
acoustic signals is a promising method to obtain the time 
- frequency energy distribution of a signal. 
These applications are based on the evaluation of the 
frequency dependent arrival time of a pulse in dispersive 
media. The underlying concept of this method will be 
briefly explained for a one-dimensional structure (e.g. a 
beam). 

A fundamental difference between most waveforms in 
structures and fluids is the dispersion. A pulse propagat- 
ing in a structure with the frequency dependent group 
velocity c g changes its shape. Due to this dispersion the 
pulse is not recognizable with correlation techniques that 
can be useful for locating airborne sound sources^. 
The wavelet transform is very useful to extract exactly 
the arrival time t a of an pulse in a dispersive media 



t a X I Cg 



(1) 



The continuous wavelet transform of a function y is 

oo 

1 



y(t) ^ [ —) dt . ( 2 ) 
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The analogue to the Fourier transforms spectrogram is 
the scalogram defined as |W^(a,6)| 2 . It can be shown 
that for a fixed scaling parameter a the arrival time t a is 
the point in time where the maximum of the scalogram 
occurs 



■>\-\Wl{a,b = t a )\\ (3) 



max (\W*(a,b)\ 

To locate a source one needs 

1. the point in time the pulse occurred, the group ve- 
locity and a sensor, or 

2. two sensors, or 

3. one sensor measuring two distinguishable wave 



types 



,14 



If the position of the source is known, it is possible to 
extract material parameters^ -. 

To improve this method, dispersion based transforms 
have been propose d 15 ' 16 , which is based on a method 
called Chirplet transfor m 17 ' 18 . These transforms im- 
proved the analysis. Nevertheless the bending wavelet 
that is presented is a new approach. 
Here a different wavelet suitable for bending waves which 
can be modelled with the Euler beam theory is presented. 
The underlying concept is not to measure the arrival time 
but to extract directly the dispersion of the pulse. The 
dispersion of the pulse is dependent on the distance be- 
tween source and receiver and the material properties. If 
it is possible to extract exactly the spreading of the pulse 
one has directly the distance or the material properties, 
depending on which is known. To define a wavelet that 
extracts the dispersion characteristics it is necessary to 
know the impulse response function in the time domain. 
For plates which can be modelled with the Euler beam 
theory this function is derived first by Boussinesqi 9 - and 
can be found in textbooks^ -. For beams it is treated in 
a companion publication- 1 since only the Green's func- 
tions for a initial deflection and velocity are found in the 
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literatur o 22 i 23 . 

The velocity v(r, t) resulting from the bending wave prop- 
agation on a infinite plate of a force impulse F a (t) = 
F S(t) at r = is 



v(r, t) 



F 



: sm 



r2\ 
4CV 



(4) 



where r is the distance from the source, C, = W B' /m", 
bending stiffness B' = Eh 3 / (12(1 - v 2 )), E the elastic or 
Young's modulus, h the plate thickness, v the Poisson's 
ratio and m" the mass per unit area. 
The bending wave velocity v(x, t) on a infinite beam, re- 
sulting of a force impulse F a (t) = F 6(t) at x = is 
given by 



v(x, t) 



F H(t) [~2 



x^_\ 
4Ci) 



(5) 



based optimisation is not reliable in finding a global op- 
timum. A second problem is the localisation of several 
overlapping pulses. A well known method that is able to 
fulfill these requirements are genetic algorithms. 



II. BENDING WAVELET 

Several different definitions based on the Morlet 
wavelet and the Chirplet transfor m 17 i 18 have been inves- 
tigated. For brevity an extensive discussion about the 
different efforts is omitted. The details of the mathemat- 
ical background of the wavelet transform can be found in 
the literatures 2 ^. 

The section begins with the definition of a wavelet with 
compact support and zero-mean. It follows a comment 
on the amplitude and frequency distribution and ends 
with possible optional definitions. 



therein x is the distance from the source, C = \J B/m', 
where B is the bending stiffness of the beam and ml mass 
per unit length. 
The term 



< 



(G) 



is the factor that controls this spreading and is called 
dispersion factor. Whereas the dispersion factor is a time 
value the nondimensional term 



Di 



% fv 



(7) 



is called dispersion number. The applicability of the 
following method depends on the dispersion number. 
Higher dispersion numbers result in a longer useful time 
period and this case is better to analyse. An exact quan- 
tification is given in the following, equation (|25|) . A high 
dispersion number is the reason for choosing a thin plate 
and a slender beam. 

In the following a new adapted wavelet will be derived 
to extract the dispersion factor from the measured pulse. 
Usually a wavelet is designed to localize a certain fre- 
quency. In contrary the proposed wavelet localizes a fre- 
quency range that is distributed over the wavelet length 
just like equation (U) or Such a choice follows the 
paradigm of signal processing, that "the analysing func- 
tion should look like the signal" . 

One may interpret the continuous Wavelet transform as 
a cross-correlation of y and ip. Hence, the idea is to find 
the function which is highly correlated with the impulse 
response. The difference is the role of the scaling param- 
eter a. It is vital to produce the presented results to use 
the scaling parameter as it is defined in equation 
The dispersion factor is determined by the scaling factor 
with the highest value of the scalogram. In principle this 
can be done with a fine grid of (a, b) values. A more ef- 
ficient way is to use an optimisation scheme. Gradient 



A. Definition 



The mother wavelet 



for train ^ 

otherwise 



(8) 



is called bending wavelet. 

A wavelet %jj must fulfill the admissibility condition 
< cty = 2tt / du> < ao. 



(9) 



where is the Fourier transform of the wavelet. The 
proposed wavelet |(5J) has a compact support (t m i n , t max ), 
which means that the admissibility condition is fulfilled 
if 



ip(t) dt = 



(10) 



holds. To fulfill the admissibility condition t m i n and t ma x 
are defined so, that equation (fPU)) holds. With the inte- 
gral - sine function Si(a;) = J sm(t)/tdt one finds that 



sin(l/i) 



dt = Si 



Si 



(11) 



Since limj^o Si(l/t) = tt/2 and that the Si-function for 
t < 2/ ir oscillates around 7r/2, one is able to chose t min 



Si(l/t„ 



S\{l/t max ) = tt/2. 



(12) 



This is a very easy option to define a wavelet and it will 
be used later to define similar wavelets. Equation (fT2"l) 
can only be solved numerically so a good approximation 
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should be used that leads to a simple expression for c^. 
The support of the wavelet is defined by 



1 



(^Tlrnax / min ljTT 



t, 



min/max 



(13) 



The function Si(l/t) and proposed possible values of t m i n 
or t m ax are plotted in figure [1] In the worst case for 
tmax — 2/7T and t m i n — > oo the difference in equation (fTU|) 
is around 0.2, but for higher values of n m i n the magnitude 
is in the order of other inaccuracies, so that it should 
be negligible. Like the Morlet wavelet, which fulfills the 
admissibility condition in a asymptotic sense the bending 
wavelet fulfills the admissibility for limi — * 0. The value 




FIG. 1. Si(l/t) and circles at l/t = (2n - 1)tt/2 



of the constant is calculated c^p with the norm in the 
Lebesque space L 2 of square integrable functions 



1/2 



\\m\h = / m 2 dt 



(14) 



The integral in equation (fl"4"|) is 
tm fsin(l/t) 2 /t 2 dt = 



4 



(15) 



With the proposed choice of t m i n and t max , the sine van- 
ishes and a normalised 1 1 -0 (i ) 1 1 2 = 1 wavelet is obtained 
if c^p is chosen to be 

__ 1 / 1 1 \ _ 7T 

* \ tmin tmax 1 ^ 



B. Displacement-invariant definition 

Wavelets that are defined by real functions have the 
property that the scalogram depends on the phase of 



the analysed function. Wavelets that are complex func- 
tions like e.g. the Morlet wavelet are called displacement- 
invariant. A wavelet tp = ip c + iip s that consists of a sine, 
tp s equation ([5]), and a cosine wavelet which is 



s(l/t) 



for tfYiiri t <^ t Ti 

otherwise 



(17) 



can be beneficial. With the integral - cosine function 
Ci(ir) = cos(t)/t dt one finds that 



Ci 



(18) 



The analogous definition of the value ir/2 for the Ci- 
function is 



Ci(l/i mm ) = Ci(l/t r , 
The approximation is given by 



0. 



1 



^minf max 



(19) 



(20) 



The effect is that the real- and the imaginary part of the 
resulting wavelet do not share the same support. This 
is an awkward definition of a mother wavelet but the 
difference between the two supports is rather small if the 
same value for n max / min is used. To keep things simple 
only the real valued sine wavelet is used in following. 



C. Orthogonality of the bending wavelet 

The trigonometric functions that are used for the 
Fourier transform establish an orthogonal base. Hence, 
the Fourier transform has the convenient characteristic 
that only one value represents one frequency in the anal- 
ysed signal. Every deviation of this is due to the window- 
ing function that is analysed with the signal. Already the 
short time Fourier transform is not orthogonal, if the dif- 
ferent windows overlap each other. Because of this over- 
lap the continuous wavelet transform can not be orthog- 
onal. The proposed wavelet should still be investigated 
since it is instructive for the interpretation of the results. 
The condition for an orthogonal basis in Lebesgue space 
L 2 of square integrable functions is 



Ojk- 



(21) 



Two different wavelets ipj and ijjk can be obtained by 
using different scaling parameters a and/or different dis- 
placement parameters b. Here the effect of two scaling 
parameters is investigated, so the following integral is to 
be solved 



ajdk sin(aj/£) sin(afc/i)/i 2 dt. 



(22) 
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To illustrate the integral, two different versions of the 
wavelet are plotted in figure [2] 



III. CONTINUOUS WAVELET TRANSFORM WITH THE 
BENDING WAVELET 
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FIG. 2. Bending wavelet (|8]) for cii = 1 (solid) and £12 = 2 
(dashed), n min = 4 and n max = 12 



One finds that 

oo 

J a,jak sin(oj/i) sin(ak/t)/t 2 dt = 

— oo 

ffi (a, eos ( a -t) sin (f) - a k cos (f) sin (*) ) 

The sine term in equation (|23p vanishes since t min / max , 
also scale with a, but actually there are two different 
values of a and so not all four sine terms vanish. With 
this result one expects a rather broad area in (a, b) with 
high values of the scalogram. 

D. Time amplitude/frequency distribution 



The application in the given context is to extract pre- 
cisely the scaling factor with the highest value. How this 
is achieved will be discussed in the next section. Here the 
realisation of a transform with a set of scaling factors is 
presented since it is illustrative. 

The algorithm implementing the continuous wavelet 
transform with the bending wavelet can not be the same 
as the algorithm implementing a transform with any con- 
tinuous wavelet, like the Morlet wavelet. The bending 
wavelet has a compact support, which must be defined 
prior to the transform. This can be done with a estima- 
tion of the frequency range and the dispersion number. 
With the equations (JT3J) and {25) it holds that 



= floor 



ceil 




(25) 



where floor (•) rounds down towards the nearest integer 
and ceil(-) rounds up. The knowledge of a useful fre- 
quency range should not provide any problems. But to 
have to know beforehand which dispersion number will 
dominate the result is rather unsatisfactory. A more 
practical solution is to calculate the corresponding n- 
value within the algorithm, which is an easy task since 
a = di. The problem with this possibility is that the sup- 
port of the wavelet changes within the transform. Since 
the support is part of the wavelet this means that strictly 
one compares the results of two different wavelets. Since 
the wavelet is normalised the effect is rather small, but 
nevertheless it should be interpreted with care. 



The time frequency distribution of the wavelet for a 
certain scaling factor a is determined by the argument 
of the sine function. The instaneous frequency ui{t) can 
theoretical be obtained with a relationship for almost pe- 
riodic functions, see Bochnei— . The actual frequency of 
the wavelet is given by 

sin p(t)-w(t) = (24) 

The \jt leading term affects the amplitude of the wavelet. 
Usually it is desired that the whole signal contributes 
linearly to the transform. To achieve this it is useful to 
have an amplitude distribution over time of the wavelet 
that is reciprocal to the amplitude distribution of the 
analysed signal. Since the bending wavelet has the same 
amplitude distribution as equation (TjQ) this may lead to 
stronger weighting of the early high frequency compo- 
nents of the impulse response. 

A force impulse that compensates the amplitude distri- 
bution of the impulse response follows a ~ 1/lo 2 depen- 
dence. 



A. Example 

To illustrate the use of the proposed wavelet the fol- 
lowing function is transformed 

/a J £sin(<z/£), for t m in <C t <C t max (ocV\ 

' ~ y otherwise, ^ ' 

with a = 10. The sampling frequency is f s = 2 9 , t m in 
and t max are defined by the corresponding values of 
fmax = fs/S and f min = 4, for convenience the point 
t = tmin is shifted to lOAi. The example function in 
equation (j2l))) is plotted in figure [31 
The example function is transformed with the algorithm 
that calculates n m i n and n max with the corresponding 
value of f m ax, fmin and a. The choice of the frequency 
range is critical. If it is too small, information will be 
lost and if it is too big the parts that overlap the pulse 
may distort the result. Here the same frequency range 
as the analysed function is used. 

The resulting scalogram is not plotted directly against 
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FIG. 3. Analysed example function (|26p 



the factor b, but shifted with the value of t m i n . This 
means that the maximum value is at 10 At (figure [4]), 
which is the value of t where f max is located. The 
maximum value is shifted when f max = / s /12 is used, as 
can be seen in figure [5] 




t/At 

FIG. 4. Contour plot of the scalogram build with the bending 
wavelet transform of equation p6[) with fmax ~ fs/8 and 

fmin — 4 

One may recognise that there are very high values if 
the wavelet is shifted and scaled along the curve a/t. This 
is expected theoretically, as discussed in section ULCl and 
can be interpreted descriptive since the wavelet does not 
localize one frequency, but has a wide frequency range 
that spreads over time. It can be quantified with equa- 
tion (|2ip . Evaluating this integral numerically for the 
values of a\ — 10, 02 = 11.75 and 62 = 7 At results in 
value of 0.68, which means that the peak at 02 = 11.75 
has 68% of the peak at a\ = 10. 

This problem of non-orthogonality is addressed by the 
following algorithm. The pulse is extracted from the sig- 
nal by first locating the position t start in the signal, where 



15 



a 10 




?10 115 120 125 130 

t/At 



FIG. 5. Contour plot of the scalogram build with the bending 
wavelet transform of equation ([26} with f max = / s /12 and 

fmin — 4 



fmax has its maximum. This is done by a Morlet wavelet 
transform with which one may find the value of t s tart 
that has the highest value of fmax- Now the transforma- 
tion with the bending wavelet is only done in the vicinity 
of t s t ar t. Technically the displacement parameters b are 
defined with t start- 

IV. LOCALIZATION WITH A GENETIC ALGORITHM 

Genetic algorithms (GA) form a particular class of evo- 
lutionary algorithms that use techniques inspired by evo- 
lutionary biology such as inheritance, mutation, selec- 
tion, and crossover. Genetic algorithms are categorized 
as global search heuristics. Details of the method can be 
found in the extensive literature, e.gi^l. 
The genetic algorithm is chosen, since it is usually very 
reliable in finding a global optimum and its ability to find 
Pareto optimums to locate several pulses. However the 
drawback is the slow convergence, that can be improved 
with a local search method. Recent publications on this 
topic ar o 28 i 29 . 

The implementation is done with functions provided by 
the open source Matlab toolbox^, if not stated other- 
wise. Principally possible but not used in the exam- 
ple is the localization of two pulses that are overlap- 
ping. For the sake of brevity a discussion on how this 
can be achieved will be omitted. The algorithm works 
with two variables, the displacement parameter and the 
scaling factor of the bending wavelet. 
The displacement parameter is defined by t star t i A/2 
or smaller values. This depends on the size of the Mor- 
let wavelet. For discrete functions it is an integer value 
but nevertheless implemented as a floating point number, 
because of lacking support for such a combination in the 
toolbox. This fact is taken into account when calculating 
the wavelet transform. A pseudo-code that describes the 
genetic algorithm can be found in the appendix [A"l 
In the end only the fittest individual is extracted. The 
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algorithm is usually quite reliable. Since it is a stochastic 
method, it can be beneficial to restart the whole process 
or to work with several sub-populations. 



A. Example 

As an example, the already transformed equation (|26p 
is investigated. The frequency range is the same as the 
example plotted in figure H As a first step the value 
tstart is calculated with a Morlet wavelet transform at 
/ s /12. The result is plotted in figure El where the max- 
imum is at U = 122. From figure [5] one may conclude 
that the correct value is 120 this slight deviation is due 
to the fact that the Morlet wavelet has a rather broad 
frequency resolution and the amplitude of the signal is 
increasing with time. If the signal sin(a/i)/i is used, the 
maximum is located at 118. The number of individuals 
is chosen to 100 and the number of generations to 400. 
The initial chromosome has a time index range of 
(106, 138) which is 3/4A and the range of the scaling fac- 
tor is (5,15). The obtained scaling factor is a = 9.99. 
The frequency range of the bending wavelet is shorter 
than the values f wmax — 37.7 < 42.7 = / s /12 and 
fwmin = 4.756 < 4 < fmin, this due to equation [25] 
The wavelet with the best scaling factor and the exam- 
ple function are plotted together in figure [7l One may 
recognize that the frequency-time distribution of both 
function match. 

3i . 1 1 . 1 



g 1.5 
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FIG. 6. Morlet wavelet transform of equation (|26[l with / = 



B. Experimental results 

Experiments on a beam and a plate are discussed 
in detail in combination with the impulse response 
functions^. The method was used to extract the disper- 
sion factor di . The theoretical functions (0} and §5§ with 
the obtained dispersion factor di are compared with the 




200 300 400 500 600 

t in [ms] 

FIG. 7. Bending wavelet at a — 9.99 and the example function 



measured curves. A good agreement of theory and exper- 
iment shows the applicability of the presented transform. 



C. Comparision with Morlet Wavelet 

The drawback of the proposed method is that it is only 
applicable to bending waves that can be modelled with 
the simple Euler beam theory or any wave that has a 
dispersion relation that follows a ~ l/V^ dependence. 
Another precondition is the rather high dispersion num- 
ber. Nevertheless bending waves are dominant in struc- 
ture borne sound problems and for thin structures the 
simplifications of the Euler bending theory are usually 
valid. 

The method is accurate, fast and easy to implement. In 
the experiments a source could be localized with a devi- 
ation lower than 10%. One can build real time applica- 
tions for source detection. The method has also principle 
advantages, it is possible to 

1. obtain the distance of a impulse or the material 
properties with only one measurement, 

2. analyse two overlapping pulse, which is not possible 
with the maximum of the Morlet wavelet transform 
where for each frequency one maximum value is 
extracted. 



V. CONCLUDING REMARKS 

A definition of a new adapted wavelet, the bending 
wavelet, is given. The mathematical properties of the 
bending wavelet are discussed. It is shown in examples 
that the transform is useful to analyse a flexural impulse 
response. 

Besides source localisation a possible application is the 
measurement of material properties in a built-in situation 
of finite structures, since the method does not depend on 
the boundary conditions. 
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The choice of a useful frequency range can be problem- 
atic. It may be useful to first analyse the signal with a 
Morlet wavelet transform to find a useful range. 



APPENDIX A: APPENDIX 

1 Linear distributed 

initial chromosome of t and a. 

2 Bending wavelet transform with 
the initial chromosome. 

3 Assignment to the 
current population. 

Repeat 

4 Self-written rank based 

fitness assignment (current pop.). 

5 Selection with 

stochastic universal sampling. 

6 Recombination with the 
extended intermediate function. 

7 Real-value mutation with 
breeder genetic algorithm. 

8 Bending wavelet transform with the 
selected chromosome. 

9 Self-written fitness based insertion 
with 70°/ new individuals. 

Until max number of generations 
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